Stability of the homogeneous Bose-Einstein condensate at large gas parameter 
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The properties of the uniform Bose gas is studied within the optimized variational perturbation 
theory (Gaussian approximation) in a self-consistent way. ft is shown that the atomic EEC with 
a repulsive interaction becomes unstable when the gas parameter 7 = pa^ exceeds a critical value 
7crit ~ 0.01. The quantum corrections beyond the Bogoliubov-Popov approximation to the energy 
density, chemical potential and pressure in powers of ^ expansions are presented. 
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I. INTRODUCTION 

The long wait after it's prediction (more than 70 years) 
for realization of the Bose-Einstein condensate (BEG) is 
possibly related to the meta-stability of the initial Bose 
gas. In fact, we require first, an atomic system would 
stay gaseous and meta-stable at very low temperature 
all the way to the BEG transition, and secondly, devel- 
opment of cooling and trapping techniques to reach the 
required regimes of temperature and density . Glearly, 
without a proper cooling technique, any ordinary atomic 
gas would undergo into a liquid or a solid state at low 
temperatures, so a meta-stale state could be created only 
with low pressure and weak interaction between atoms. 

Even once created, the condensate still remains as a 
fragile and subtle object 0. The enemies of BEG such 
as crystallization, disassociation and three-body recom- 
bination may easily destroy it within very short time. 
When the sign of interaction (or equivalently of the s- 
wave scattering length, a ) is suddenly changed into a 
negative value, the BEG collapses and then undergoes 
an explosion in which a substantial fraction of the atoms 
were blown off [Bosenova) 3, 4]. 

Due to the " bad collisions" , even an atomic BEG with 
a repulsive interaction has a limited life time. Recently, 
Cornish et. al. 5] carried out an ingenious experiment 
with spin polarized atomic ®^Rb. In the experiment, they 
showed that one could control the strength of interatomic 
interaction for the BEG by employing the Feshbach reso- 
nance method. A very large value of the scattering length 
(a « 4500 A) has been achieved in this experiment, which 
corresponds to the gas parameter of the condensate to be 
about ^max ~ 0.01. This phenomenon has been recently 
studied by Yin in Random Phase Approximation. The 
author has shown that when 7 exceeds a certain criti- 
cal value the molecular excitation energy becomes imagi- 
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nary and, hence, the atomic BEG is dynamically unstable 
against molecular formation. 

It is well known that many of the basic properties of 
the condensate of dilute Bose gases in existing experi- 
ments can be described reasonably well using the mean 
field approximation (MFA) which reduces the problem to 
the classical Gross - Pitaevskii equation (GPE) f?\. How- 
ever, fluctuations of the quantum field around the mean 
field provide corrections which become increasingly im- 
portant as higher condensate densities (say large 7 ) are 
achieved. It is therefore important to understand the ef- 
fects of quantum field fluctuations especially at large gas 
parameters. 

In the present paper, we study the properties of a ho- 
mogeneous atomic Bose gas using optimized Gaussian 
approximation It has been proven that the corre- 
sponding Gaussian effective potential contains one loop, 
sum of all daisy and superdaisy graphs of perturbation 
theory [Qj] and leading order in 1/A'^ expansion. 

The first application of the Gaussian variational ap- 
proach to a uniform BEG was done by Bijlsma and Stoof 
ten years ago [l3|. However, it was pointed out in ex- 
cellent review by Andersen that flTl even a modified (by 
introducing many body T-matrix) Gaussian approxima- 
tion of Ref. [l3| does not satisfy the Hugenholtz-Pines 
(H-P) theorem especially at very low temperatures. This 
is particularly caused by a long standing problem encoun- 
tered in the most of field theoretical approximations: it 
is impossible to satisfy the H-P theorem, namely mak- 
ing the theory gapless at the same time and maintaining 
the number of particles with the same value of the chem- 
ical potential. In other words, the chemical potential 
defined by the H-P theorem does not coincide with the 
chemical potential found from the minimization of the 
thermodynamic potential with respect to the condensate 
density. Note that, even the T-matrix approximation 
cannot resolve this problem completely since in this case 
one gets "mismatch of approximations" which makes the 
approach as non self-consistent. 

One of the possible solutions of the above mentioned 
problem has been proposed recently by Yukalov 12]. He 
has shown that Hartree-Fock approximation (HFA) can 
be made both conserving and gapless by taking into ac- 
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count of two normalization conditions instead of one. 
Hence, two chemical potentials each for the condensed 
fraction (^o) and the uncondensed fractions (/ii) should 
be introduced to describe the BEC self-consistently. 

In the present paper, we reformulate the field theoret- 
ical Gaussian approximation following the Yukalov pre- 
scription, and apply this self-consistent approach to in- 
vestigate the properties of a uniform BEC. 

The paper is organized as follows: In Sect. II, we 
extend the field theoretical approach by implementing 
Yukalov's ideas. In Sect. Ill, we calculate the free energy 
in a Gaussian approximation, and also show it's relation 
to the one loop and Bogoliubov-Popov approximations 
(BPA). In the next two sections, we present the proce- 
dure of minimization of the free energy. The numerical 
results and their discussion are presented in section VI. 
Sect. VII summarizes the paper. 



II. QUANTUM FIELD FORMULATION WITH 
YUKALOV PRESCRIPTION 

A grand canonical ensemble of Bose particles with a 
short range s - wave interaction is governed by the Eu- 
chdian action [ui, [131 

S[<P,r]= f^dr j dr{r{T,r)[dr - -M]V'(T,r) 

-H|r(r,r)V'(r,r)]n,(l) 

where V'* {T^ f) is a complex field operator that creates a 
boson at the position r, /i the chemical potential, g the 
coupling constant given by 47ra/m, m the atomic mass 
and (3 = 1/T the inverse of temperature T. The free 
energy of the system can be determined as 



= -T\nZ, 
where Z is the functional integral, 

VipV'^' e-xp{-S[ip,i)*]}, 



(2) 



(3) 



performed all over Bose fields tp and tp* periodic in 
T g [0, /?]. When the temperature in a Bose system falls 
below the condensation temperature T^, breaking of the 
?7(1) gauge symmetry may be taken into account by the 
Bogoluibov shift of the field operator. 



i>{r, r) — v{t, r) + tp{T, r), 



(4) 



where u(t, r) is the condensate order parameter. In the 
uniform system v{T,r) is a real constant, v{T,r) — v, 
ip{T, r) is the field operator of the uncondensed parti- 
cles satisfying the same Bose commutation relation as 
ip(T,r). The conservation of particle numbers requires 
that tpir, r) has non - zero momentum component so that 



and ip and v are orthogonal each other 
dfili{r)v{r) = 0. 



(6) 



The condensate order parameter v defines the density of 
condensed particles while ip defines the density of uncon- 
densed particles: 



po , pi =< tjj* {r)'ip{r) > 



(7) 



Having performed the Bogoliubov shift, one may in- 
troduce the grand canonical thermodynamic potential of 
the system as 

n{p,v) = T{p,v)\^^^^y (8) 

In a stable equilibrium, Q attains the minimum: 



dn{p,v) d^n{p,v) ^ ^ 



dv 



d^v 



(9) 



Apart from the H-P theorem, the chemical potential 
should satisfy the normalization condition 



N = { drij}* {r)'4){r)) , 



(10) 



where N is the total number of particles. However, as 
it was pointed out in the above, the chemical potential 
corresponding to the minimum of may not correspond 
to the chemical potential p, determined from the normal- 
ization condition. 

To overcome this difficulties, Yukalov proposed to 

• Introduce one more normalization condition A^q = 
PqV . So that for the uniform system 

No+Ni=N, Ni^iJ dr^*{r)^{r)), (11) 

which simply states that the total number of par- 
ticles should be equal to the sum of the number of 
condensed and uncondensed particles. 

• Introduce two chemical potentials /ig and pi^ for 

the condensed and the uncondensed fractions re- 
spectively as well as a Lagrange multiplier A to 
satisfy the Eq. The total system chemical po- 
tential, p — —{dfl/dN), is given by 



N 



(12) 



These prescriptions lead to the following action. 



S[^,r] = I drj df{riT,r)[dr - ^]V'(T,r) 
-Pi'P*{t, r)ip{T, r) - pov^ - Aip{T, r) 



-A*r{r,r) + ^[r{T,rMT,r)]'}, (13) 

which should be used in Eq. ([3]). Further, po can be 
determined from the minimum condition Eq. ([9|), while 
Pi itself by the requirement of H-P theorem 



W = 0, 



(5) 



Pi — Sii — El 



(14) 
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where Sn and S12 are the normal and the anomalous 
self-energies. As to the condensed fraction TVq, it could 
be found by solving the normalization Eq. (jlip where 
the uncondensed fraction TVi is given by 



iVi 



i9/ii 



(15) 



III. GAUSSIAN, ONE - LOOP AND 
BOGOLIUBOV-POPOV APPROXIMATIONS 



In the present section, we show how this scheme can be 
realized in practice. Substituting Eq. (j4|) into Eq. ([TS]) . 
one may rewrite the action in powers of v and ip, 



S 
Sio) 

5(3) 

5(4) 



J d,T J dr{-u^^o + — }, 



dr J dr[gv^-A*-A]['ip*+-ip], 
dr / dr{ip*[dr - 



2 m 



f ^ 

dr / d'rtp*tp*^tp. 



(16) 



In the following, S^^'^ will be omitted since it can be set 
to zero by an appropriate choice of A in order to satisfy 
Eq. (0). 

Now in accordance with the variational perturbation 
theory, we add and subtract the following term: 



parts: 

S = 

c 

'-^clas ~ 

c _ 

*-^free ~ 

"^int ^ 



"^clas "^frcc '^int' 



gv 



dr / df[ieab-iJ-'adTiJb + i>ii-7, 1" 



2 m 



S. 



int 



2m 

mt mt mt 



dT / df[ijl{3gv^^Ui,) + tl4{gv^ 



(19) 



n 



22) 



mt 



:.(4) 
'int 



V2J0 



I dr I drvMi^l+4^1), 



dT / dr{ri + 'irir2 + r2)- 



(20) 



Here, tabiO'ib = 1,2) is the antisymmetric tensor in 
two dimensions with ei2 — 1 and following notations are 
introduced, 



(21) 



llii — Ell + E12, 1122 = Ell — E12, 

Xi = IIii — ^il, X2 = 1122 — Ml- 

In accordance with Refs. [HI, [HI, Ilab are the compo- 
nents of the 2x2 self-energy matrix. 

The free part of the action, S'fj.gg in Eq. ((20)) gives 
rise to a propagator, which can be used in perturbative 
framework. In a momentum space, 



V'a(T,r) 



^V'ti('-^n,fc)exp{jLj„r-|-ifef}, (22) 



where dk/{2nf, and ojn — "iimT is the Mat- 

k 

subara frequency. The propagator is given by 
1 



G(w„, k) 



El 



£k + X2 OJn 
—bJn £k + Xi 



(23) 



df 



(17) 



assuming Sn and S12 as real constants. Further, we 
write the quantum fluctuating field in terms of two 
real fields 



V> = -^(i/ji + ii/)2), it* = -^{ipi - i^p2) ■ (18) 



After some algebraic manipulations [ij, ll5[, one can 
split the action into "classical" , "free" , and "interaction" 



with the dispersion relation, = [sk + Xi){ek + X2) 
and Ek — K /2m. 

With this Green's function using Eqs. ([2]) and ([3]), 
and neglecting terms S*!^! and >S'j'j^|_, one may get the 
thermodynamic potential in the one loop approximation: 

n(i^)(^.o,Mi,»^) = ^(-Mo^' + ^) 



(24) 



k k 

+ i [B(35«' - nil) + A{gv^ - n22)] 



with IIii = "igv^ and 1122 = QV^ ( ^ and B will be given 
below), so that the last term in square bracket can be 
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dropped. Note that, hereafter we perform explicit sum- 
mation by Matsubara frequencies (see e.g. [lal)- As to 
the BPA, it can be obtained by introducing an auxihary 
expansion parameter as it was shown by Kleinert 

Loop expansion of may be organized by using the 
propagator G{ujn,k) with constraints Xi = 2gv^ and 
X2 = as illustrated in Ref. To take into account 

higher order quantum fluctuations, one has to calculate 



^"^int^ and Although these quantities can not be 

evaluated exactly, they may be estimated in the Gaussian 
approximation [25j . where for the homogeneous system: 



IV. THE GAP EQUATIONS AND THE 
THERMODYNAMIC POTENTIAL AT T = 



In this section, the variational parameters IIn and 1122 
will be determined using the principle of minimal sensi- 
tivity 3. From Eqs. ((^ and (P7)l . the gap equations 
may be found 

dQ{Xi,X2,v,iJ.o,Hi) 
dXi 

= l{A[[9v^ - Ml - X2] + Bi[3gv^ - Ml - Xi] (^0) 
+g[A[{3A + B) + B[{3B + A)]/2V} = 0, 



0, 



(25) 



int' 

{^UD = i^t) = 3gL(o), 

n — — 00 k 

G22(0) = V-^A. 



Finally, combining Eqs. (|24p and (|25p. we get the follow- 
ing expressions for the thermodynamic potential: 



and 



k k 

+ ^ [B{3gv^ - Uii) + A{gv^ - n22) 



(26) 



2AB] , 



where 



E 



^-E 



Ek+Xi 



1 



1 

2 + ixp(/3Sfe) - 



(27) 



£k + X2 
Ek 



2 exp(/3Sfc)-l 



The free energy in Eq. (|26p is supposed to have all the 
information about the system. Particularly taking it's 
derivative with respect to /ii, one gets the expression for 
the uncondensed fraction A'^i: 



do. 



^^{A + B- {3gv^ - Uii)B' - {gv^ - n22)A' (28) 
-^[{2,A + B)A' + {W + A)B']}, 

where A' = ^A/^^ll and B' = ^B/^^ll. Note that the 
same expression for the uncondensed fraction could be 
obtained in an alternative way as 

pi = ^ =< >= J VtjjVi)* ex.p{ — S[tl},'4!'']}'4>*'tp 

(29) 



dn{Xi,X2,v,fJ.o,iJ,i) 

8X2 

= 2'f^2[fff^ - Ml - ^2] + B'2[igv^ - Ml - ^1] 
+g[A'2{?,A + B) + B'2{3B + A)]/2V} = 0. 



(31) 



where 



_ OA 1 _ dB _ , 

= dTi - I^W, - dT2 = 

_ dA _ l^ {ek+Xif 



as _ 1 



El 
(£fc + X2)' 



(32) 



Above, we have two equations ([50]) and (pij) with re- 
spect to three unknown quantities {Xi, X2, /iti}. An ad- 
ditional equation is supplied from the relation between 
the chemical potential and the self-energies given by the 
H-P theorem. So, from Eqs. (fT4|) and pT|) . one can 
immediately conclude X2 — 0, and, hence, in the long 
wavelength limit (fc 0), the quasiparticle energy E]^ 
behaves as ck (with c — Xi/2m ) thus being gapless, 
as expected. With this constraint, the gap equations may 
be simplified as 



2Xl-fMl-5p^;^--:^/l,l(Xl)=0, (33) 
h,i[Xi)[Zgv^ - Ml - X{\ - J_2,-i(Xi)[(7«2 - Ml] 



+ ^^'l^'\ ^hAXi)+I-2,-i{Xi)] = 0, 
and Eq. as 

iVi = i/i,i(Xi) 

[Jo,i(Xi) - 2I.2,-i{Xi)][9,Vgv^ - gh,i(Xi) - 8Vmi 
128V"m 

Here, the following dimcnsionlcss integral 



(34) 



(35) 



4.(Xi)=E 



(36) 



is introduced. Their explicit expressions and the relations 
between them evaluated in dimensional regularization are 
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presented in the Appendix of Ref. In particular, 

— 2B\x2=o ~ — 4j4|x2=o, 



- ^ 



(37) 



dXi 



Io,i{Xi) 



Note that Eqs. ^ and ^ include /_2,-i (-'^i), 
which is infrared divergent, /_2,-i(Xi) ^ 1/e + 
\nK^ /mXi (with e 0). Below we show that this in- 
tegral will be canceled exactly. In fact, eliminating Xi 
from Eq. dnSl) as 



(38) 



and substituting it into Eq. ([34|) . one observes that the 
latter is factorizcd: 

[Io,^{X^)~2I-2,-l{Xl)]lgh,liXl)-8Vgv^+8V^il] = 0. (39) 

Finally, from the last two equations, we find the formal 
solutions of the gap equations 



Ml , 5gv' llgJi,i(Xi) 



X^^2gv^ + ^h,^{X^), 



^^1 =g«'-^/i,i(Xi) 



(40) 



(41) 



We denote these optimum values of Xi and /ii by Xi 
and fii, respectively, which are explicitly dependent on 
w^. Now, comparing Eqs. ([55]) and ([55]) . one can easily 
see that only the first term in Eq. ((35l) survives 



iVi = g/i,i(Xi). 



(42) 



Now, inserting these formal solutions into Eq. ((26|) and 
using the relations between the integrals, Eq. ([ST]) , gives 
the following form for fl 



nglUXi) 



(43) 



128F 



where 



/o,-i(^i) = — > ^ VeTVEfe +Xi = J- ' . 

k 

(44) 

In particular, neglecting in Eq. the last term gives 
the one-loop result: 



n{Xl,V,fJ.o)\xi=23v2 = r2'^^>(M0,Ml,w)lMi=9i.2- 



(45) 



presented in the previous section. In the stable equilib- 
rium, the grand canonical potential reaches the global 
minimum as a function of v : 



dQ{Xi,v, fip) 



Vp{-HQ + gpna) 



dno 

4 + l6Vm J - ° 



(46) 



where no = v'^/p = Nq/N and X[ = (dXi/dno). Note 
that the same equation could be obtained from the orig- 
inal equation (^5]) as: 



dfl{Xi,X2,v, po, fJ,i) _ dQ 
dno dno 

on \ dXi 



\ djii J dno \ dXi J dno \ 8X2 



on \ 8X2 



8no 



0, 



(47) 



where the last two terms may be omitted due to the gap 
Eqs. ([30|) and (|3ip. and the factor in the second term is 
related to iVi bv p5)) . 

Clearly, the optimal value of w^, i.e. defined by Eq. 
(|46p . should correspond to the normalization condition 
in Eq. pT|) (constraint): 



-2 , /V \ -2 , Il,l{Xi) 
V +pi{Xx)=V + 



(48) 



which may be considered as a nonlinear equation with 
respect to the c- number v"^ with a fixed p and Xi{v'^). 

Strictly speaking, must be determined from Eq. 
(I46p as a function of /io, and after substituting it into 
Eq. the latter should be solved with respect to /zq. 

However, this would be a rather complicated way, since 
Eq. (Uni) is a highly nonlinear equation. On the other 
hand, one may assume that is known as a solution of 
Eq. and could be extracted from Eq. 

Following this strategy, we obtain 



Mo = gpno + 



XUi,i(Xi) 

4\//9 



llSr/o.i(Xi) 



l6Vm 



(49) 



in particular, neglecting the second term in square brack- 
ets and taking into account Xl'" = 2gpno, we have po for 
the one-loop approximation 



IL 

Mo = gp 



no + 



2Vp 



(50) 



Further simplification, by introdu cing an auxiliary ex- 
pansion parameter rjiL as in Ref. [l7|, gives po for the 
BPA 



BP 

Mo 



gp 



1 



gIo,i{Xi = 2gp) 



2Vm 



(51) 



As to the total system chemical potential /i, it follows 
from Eqs. and dH]) as 



(52) 



p = pini + pono, 



where ni = 1 — no, and pi and po are given by Eqs. (HT|) 
and (149]) . respectively. Now, substituting (|49|) into (|43|) . 
one may obtain the pressure as P = —fl/V 



P = |ffngp2 + _[noX;/i,i(^i) - 2m/o,-i(Xi) 



+ "l28ml^'^ [2noX{/o,i(Xi) +mJi,i(XO 



(53) 



The ground state energy density of the BEC, £, may 
be obtained by a well known formula £ = {il + pN)/V. 
This may be easily done by rewriting the term pqv'^V in 
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Eq. pS)) as ^qv'^V = fj,N — fiiUiN (which follows from 
Eq. dnni)) and using Eq. dU): 



lo 2 2 (54) 

+ ^Ji,i(Xi)-^/i,i(XO. 



Now, one may immediately obtain 

+ + 12SV 



and 



^ - — + 2Y^0,-l(^l)\Xi^2gp, 



(55) 



(56) 



for the Gaussian and Bogoliubov - Popov approximations 
respectively. 

It is well known that in the BPA, the normal (Sii) and 
the anomalous self-energies (S12) are rather simple (l7j : 



= 2p/5, Sfa^ = gp. 



(57) 



In the Gaussian approximation using Eqs. (PT|) . PO)) . 
and (1411). one obtains 



2 



Si2 = :ii = <;.V|^/M(Xi) 



4V 



(58) 



which can be further simplified at the stationary point as 

Sii = 2gp 
E12 =3p(l + 2ni). 



(59) 



Clearly, neglecting the uncondensed fraction ni in the 
last equation, we recover the Bogoliubov - Popov approx- 
imation, (|57p . The dimensionless sound velocity defined 



as c = lim Ek/k ~ y Xi/2m is simply related to S12 as 



(60) 



V. SOLUTIONS TO THE GAP EQUATIONS 

In this section, we analysis possible solutions to the 
gap equation which can be written as 



(61) 



Before solving this equation, we emphasize that in accor- 
dance with the general principle of the variational Gaus- 
sian approximation, the constraint in Eq. (|48p and the 
procedure of minimization of the free energy with respect 
to v'^ may be imposed only after finding an explicit ex- 
pression for Xi = Xi{v'^) as a function of w^, which can 
be done by solving Eq. ([6T|) analytically. Note that when 
the second term on the RHS of Eq. is neglected, one 



obtains a well known result of the one-loop approxima- 
tion: Xl^ = 2gv'^, and further, assuming here = p 
gives the self-energy for BPA : X(^^ = 2pg. 

In general, the Eq. (j6ip can be rewritten in a dimen- 
sionless form 



432^ 



3456 



3/2 



(62) 



where the following dimensionless quantities were intro- 
duced 

Z^^X,/2gp, iV, = l^, (63) 

TT 

with 7 — a^p is the gas parameter. Analysis shows that 
Eq. ([62]) has no real positive solution when > 1. This 
is illustrated in Fig. 1 where the solid curve presents 
RHS, and the dashed straight lines present LHS of Eq. 



(162)) for = 0.1; 0.3; 0.7; 1.0; 1.1 from the bottom to 
the top, respectively. It is seen that when Nj < 1, there 
are two different solutions (denoted as crosses in Fig.l 
) which overlap at ^ 1 and Z = tt/IU = 0.0218, 
and then disappear. This is one of our main results 
confirming that there is a critical value of 7, or more 
exactly critical value of Nq'j/N which controls the sta- 
bility of the uniform Bose condensate at T = 0. When 
N~f = 432no7/7r exceeds unity, {N-^ > 1), Xi and hence 
the self-energy becomes complex, and the BEG will be 
unstable. 

Differentiating Eq. ([62]) by and solving with re- 
spect to dZ/dN^, one obtains: 



Z' 



dZ 



^3/2 



dN^ 432(0F- 12\/Z)' 



(64) 



which is singular at Z — 7r/144, i.e., at = 1. Thus, at 
the critical point, — 1, 



lim — — = 00, 



(65) 



and hence, at this point the chemical potential of the 
condensate fiQ in Eq. which is responsible for the 

thermodynamical stability of the system, has a singular- 
ity- 

For N-y < 1, the solutions are given as [2y| 



Zi = 



Z-2 - 



— [2C1 C0S(C2) 



TT 

64 



216 



■3] 



(66) 



-^[-ci cos(c2) -t- VSci sin(c2) + 3] (67) 
576 



I32" 



%/37r7V. 



3/2 
7 



1944 



1944 



(68) 



where ci = ^9 - %N^, ci = arccos{[27 - Z^oN^, -|- 
87V^2]/c?}/3. 

It is understood that only the second solution, Z2, is a 
physical one, since for the case oi Z = Zi the self-energy 
Xi is irregular at 7 — s- 0. Moreover, only Z = Zi corre- 
sponds to the minimum of the thermodynamic potential. 
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FIG. 1: Graphical solution of the gap equation (|62|l . The solid 
curve represents RHS, and the dashed straight lines represent 
LHS of the equation for iV.^ = 0.1; 0.3; 0.7; 1.0; 1.1 from the 
bottom to the top, respectively. 



FIG. 2: Condensate fraction no = 1-0(7) ^ ^ function of 
7 = pa''. Solid, dashed and dot-dashed curves correspond to 
the Gaussian, the one-loop and Bogoliubov-Popov approxi- 
mations. 



{d^n/d^Xi)\^=^^ > 0. Thus, we conclude, Xi = 2gpZ/-f 
with Z — Z2. In particular, taking into account the first 
term in the expansion of Z2 in Eq. (|68p . one obtains 
Xi « 2gpno = Xl^ as expected. 



VI. RESULTS AND DISCUSSIONS 

Expansion for small 7. The starting point of our nu- 
merical calculations is the Eq. ([5T|) . which can be rewrit- 
ten as 



1 - no 



87^/7? 



0. 



Before analyzing this nonlinear equation we note that 
the majority of experiments with ultracold trapped gases 
deal with weakly interacting atoms, so that 7 is very 
small, i.e. 7 ~ 10-9 . . . 10"''. Thus, to obtain a low- 
density expansion of physical quantities, one may search 
for the solutions of Eq. ([69]) in power series of ^/j to get 



ni 



= 1 - 

8/7\i/2 



64 7 



0(7^/2). 



(70) 



Clearly, the first term corresponds to the Bogoliubov ap- 
proximation, while the others may be considered as quan- 
tum corrections to this approximation. Note also that 
the above expansion Eq. (jTO]) is exactly the same as the 
one obtained in the modified Hartree-Fock Bogoliubov 
(HFB) approximation Now, using Eqs. jH]), (|49l) . 
(153)1 , jini), dSH), §0^ and ^ we obtain the following 
low density expansions for the energy density, chemical 
potentials, self energies, the sound velocity and the pres- 
sure: 



£^gp^l + 



128^7 ^ 1287 



S^/^T on 
(71) 



_ 16^7 1287, ^ _ /I ^ ^6^7^ 1287, 

pi ^ gp{l~ — — }, S12 ~gp{l + --^ + — — }. 

(72) 



9P\. I 64^7 647, 
2 ^ 50F Stt ^' 
(73) 



which are in good agreement with BPA ^19*1. 

Critical density and exact solutions. In order 
to discuss exact solutions of the equation, ((69)) . we first 
establish the boundary for 7 which is related to the crit- 
ical value of N^^ found in the previous section. This may 
be evaluated directly by substituting N.y = 1, Z = 7r/144 
into Eq. ([69)1 . which immediately gives jcr — 57r/1296 ~ 
0.012120. It is interesting to observe that when 7 ap- 
proaches this critical value, the condensed fraction re- 
mains still large, lim no(7) = 7i'/4327cr = 3/5 = 0.6 

but the condensate as a whole become unstable. 

Fig. 2 presents the condensate fraction ^0(7) in the 
Gaussian (solid line), the one- loop (dotted line) and 
Bogoliubov-Popov approximations. It is seen that due 
to the quantum fluctuations the condensed fraction de- 
creases faster (with increasing 7) in the Gaussian approx- 
imation than in BPA. 

The chemical potential p, = fj,ono + piUi is presented in 
Fig. 3. One may observe that, in the Gaussian approxi- 
mation it varies slowly with increasing 7 , almost coincid- 
ing with that for the BPA. However, when 7 approaches 
the critical value jcrit = 0.012 it starts to increase very 
fast since in this region, when iV^ 1, Ar( in Eq. ([49]) 
becomes very large and so does fio- Bearing in mind that 
the chemical potential is the energy needed to add (or 
extract) one more particle to (or from) the system, one 
may interpret this effect as a particle number saturation 
of the condensed particles. In other words, when 7 (or 
more exactly 7V^) reaches the critical value, the num- 
ber of condensed atoms Nq cannot be further increased. 
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FIG. 3: The chemical potential in the Gaussian (solid line) 
and Bogoliubov-Popov approximations. 



since it will lead to a dynamical instability of the BEC. 
The pressure defined by eq. ([53]) is positive in the region 
7 < Jcrit, but near the critical point 7 ~ jcrit it becomes 
negative and small, as expected. Note that at this point 
the energy density, the self-energies and sound velocity 
remain finite, since corresponding expressions Eqs. ([55])- 
([50)1 do not include X[ explicitly. 

Now we consider possible origin of the instability found 
above. It is well known that ^] the BEC is an effect of 
the exchange coupling, which leads to an effective at- 
traction between atoms forcing them to accumulate in a 
single state. However, when the density (or scattering 
length) reaches a critical value this effective attraction 
makes the condensate collapse. 

Another possible reason is a three body recombination 
of condensed atoms. Although there is no explicit 3 body 
interaction in our starting Lagrangian, it was shown that 
PH . at T ^ 0, the repulsive two body interaction leads to 
a three body recombination with the rate constant a^ec oc 
and, hence, the three body recombination becomes 
very significant for a large scattering length i.e. large 
7. This seems to be one of the main reasons for the 
fact that a stable condensate with large gas parameter is 
inaccessible experimentally. When 7 exceeds the critical 
value, the atoms start to combine into molecules and the 
condensate may undergo phase transition into a solid or 
a liquid state. 



VII. SUMMARY 

In conclusion, we have developed a new Bosonic self- 
consistent variational perturbation theory, which can be 
made the starting point for systematic expansion proce- 
dure [l^. We have shown that taking into account two 
normalization conditions at the same time solves the old 
outstanding problem of Bose systems making variational 
perturbation theory both conserving and gapless. 

Studying the properties of a system of uniform Bose 
gas at zero temperature with repulsive interaction both 
analytically and numerically, we have found that in this 
system there is a dynamical parameter N-y oc no pa^ which 
controls the stability of the Bose condensate. When this 
parameter remains smaller than the critical value the 
phonon spectrum is purely real and the excitations have 
infinite lifetimes. On the contrary, when exceeds the 
critical value the condensate becomes unstable, in similar 
fashion to the BEC with an attractive interaction. Note 
that this phenomena cannot be obtained in ordinary per- 
turbative framework. 

It would be quite interesting to study the dependence 
of critical on temperature. It was discovered long 
ago by Bethe [2^ that the inelastic cross section, which 
tends to destroy the condensate, varies as 1/fc, (called 
as 1/velocity law), and hence the bad collisions can be 
surprisingly large near zero temperature. Thus, the tem- 
perature dependence of the critical seems not to be 
trivial. This work is on progress. 
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